Stability of exponentially damped oscillations under perturbations of the Mori-Chain

The preliminary section already touched on the relation between a correlation function C(t) and the corresponding Lanczos coefficients b_n . As mentioned, this correspondence is one-to-one, thus, a set of b_n 's uniquely determines C(t) and vice versa. However, the convergence can be quite subtle, i.e., there may be fairly similar dynamics with vastly different Lanczos coefficients. On the other side, similar Lanczos coefficients may lead to quite different dynamics.

In the numerical experiments below, we proceed as follows. First, we choose a correlation function $C^{\prime} (t)$ whose stability we want to probe. We may specify $C^{\prime} (t)$ as a concrete analytical function, or as a member of a class of functions, e.g., exponential decays. The exact Lanczos coefficients corresponding to $C^{\prime} (t)$ are denoted by ${b}_{n}^{{\prime} }$. These coefficients are, in general, unknown and obtaining them is quite a difficult task. Our objective is to find Lanczos coefficients b_n corresponding to a correlation function C(t) that should be practically indistinguishable from $C^{\prime} (t)$ for all intents and purposes $[C(t)\simeq C^{\prime} (t)$]. These 'approximate' coefficients b_n may be obtained in three different ways:

• (i)  
  The exact coefficients ${b}_{n}^{{\prime} }$ may be analytically known. In this case ${b}_{n}={b}_{n}^{{\prime} }$ as well as $C(t)=C^{\prime} (t)$.
• (ii)  
  The coefficients may be obtained via an educated guess. In this case, the b_n may be quite different from ${b}_{n}^{{\prime} }$ while still $C(t)\simeq C^{\prime} (t)$.
• (iii)  
  The coefficients may be 'reverse-engineered' from the correlation function $C^{\prime} (t)$. For details see appendix A.

In the following, we will omit the technical distinction between C(t) and $C^{\prime} (t)$, since we assume that Lanczos coefficients can be found for which the two correlation functions are practically indistinguishable.

In the following, we compare dynamics of correlation functions C_A (t) and C_B (t) with respect to their stability under a certain class of perturbations. Corresponding sets of Lanczos coefficients b^A _n and b^B _n may be obtained by one of the three methods mentioned in the last section. We demand that these Lanczos coefficients satisfy the following [partly physically motivated] requirements sufficiently well:

• (i)  
  The Lanczos coefficients b^X _n should reproduce the respective correlation function C_X (t) to an acceptable accuracy, where X = A, B (this is sort of obvious and has already been addressed above).
• (ii)  
  The Lanczos coefficients b^X _n should comply with the universal operator growth hypothesis, i.e., they should eventually attain linear growth.
• (iii)  
  The resulting correlation functions C_A (t) and C_B (t) should decay on more or less the same time scale. This is to ensure a fair comparison between the two, since faster dynamics are typically less affected by perturbations than slower ones.
• (iv)  
  The Lanczos coefficients b^A _n and b^B _n should be similar in magnitude. In particular, we demand that the quantity ${\sum }_{n}{\left({b}_{n}^{X}\right)}^{2}$ is comparable for X = A, B. We show in appendix B that this sum is related to the spectral variance of the Hamiltonian. Hence, this condition fosters the notion that the two correlation functions C_A (t) and C_B (t) originate from different observables while the underlying Hamiltonians are quite similar. In practice, a value of

  $$\begin{eqnarray}&&q({b}_{n}^{{\rm{A}}},{b}_{n}^{{\rm{B}}})=\displaystyle \frac{{\sum }_{n}{\left({b}_{n}^{{\rm{A}}}\right)}^{2}}{{\sum }_{n}{\left({b}_{n}^{{\rm{B}}}\right)}^{2}}\end{eqnarray} \tag{ 8 }$$

  close to unity is desirable.

The considered perturbations are designed on the level of Lanczos coefficients. The coefficients b_n corresponding to some correlation function C(t) will be slightly altered according to

$$\begin{eqnarray}&&{\tilde{b}}_{n}={b}_{n}+\lambda {v}_{n},\end{eqnarray} \tag{ 9 }$$

where the perturbed coefficients are denoted by ${\tilde{b}}_{n}$. Here, λ is the perturbation strength and v_n is specified below. We show in appendix C that this particular form of the perturbation yields a sensible scaling of the perturbed Hamiltonian with λ.

In general, it is an intricate problem to determine how a perturbation in the form of v_n corresponds to a perturbation ${ \mathcal V }$ on the level of Hermitian matrices. This is an issue that we do not attempt to tackle in this paper.

Here, we try to make as few assumptions as possible and model the perturbation v_n as random numbers [with the only restriction of a minimal correlation length, see below]. Concretely, we set

$$\begin{eqnarray}&&{v}_{n}=\sum _{k=1}^{{N}_{{\rm{f}}}}{x}_{k}\cos (2\pi {nk}/d)+{y}_{k}\sin (2\pi {nk}/d),\end{eqnarray} \tag{ 10 }$$

where the x_k , y_k are real, random numbers from a Gaussian distribution with zero mean and unit variance. They are normalized as ${\sum }_{k}{x}_{k}^{2}+{y}_{k}^{2}=1$.

The sum in equation (10) is capped at a number N_f. This corresponds to excluding shorter wavelength or higher frequencies, which induces a minimal correlation length in the coefficients ${\tilde{b}}_{n}$ [reddish noise]. No truncation at all, i.e., N_f = d, would be equivalent to adding uncorrelated random numbers [white noise] to the coefficients b_n . This choice can be justified by invoking the tight-binding model interpretation in which the Lanczos coefficients represent hopping amplitudes between neighboring sites. Simply altering the hopping amplitudes in a random, uncorrelated manner leads to localization effects similar to Anderson localization [18]. Through unsystematic varying of N_f, we find that the transition from unlocalized to to localized behavior is quite sharp. Since the aim is to study the stability of certain relaxation dynamics, we choose N_f as large as possible, but small enough to avoid said localization effects. In practice, this amounts to N_f ≈ d/3. In section 4.3, we ease this restriction to identify 'pathological' perturbations.

Once the perturbation has acted on the coefficients b_n and yielded the perturbed coefficients ${\tilde{b}}_{n}$, we can calculate the corresponding perturbed dynamics $\tilde{C}(t)$ as laid out in section 2. Now, some tools are needed in order to judge how strongly the dynamics were altered by the perturbation. In particular, it is important to assess to what extent the perturbed dynamics still falls into the original class of functions to which the unperturbed dynamics belonged. For example, we may ask if $\tilde{C}(t)$ is still exponential, given that the original dynamics C(t) was exponential [the decay constants may differ]. To this end, we attempt to fit the perturbed dynamics with a function that also described the unperturbed dynamics, e.g., we may try to fit $\tilde{C}(t)$ with f(t) = Ae^{−μ t} . In practice, the fit is obtained via a standard fitting routine.

Before we introduce quantifiers that measure the effect of the perturbation, some remarks on the fitting ansatz are in order. Correlation functions are symmetric in time, thus, all odd moments necessarily vanish. In particular, correlation functions have zero slope at t = 0. Further, we normalize the correlation function to unity at t = 0. However, for fitting we choose functions that lack these features, i.e., the above exponential ansatz may yield an A that slightly differs from unity. Further, the slope of an exponential at t = 0 is never zero. This 'negligence' is due to the fact that we are more interested in an overall description of relaxation dynamics, rather than in the details of the short-time behavior.

We assess the quality of a fit f(t) by calculating 'how far off' it is from the given perturbed dynamics. Concretely, we define a measure of the 'error' or 'deviation' by the expression

$$\begin{eqnarray}&&\epsilon =\sqrt{\displaystyle \frac{1}{{N}_{\mathrm{eq}}}\sum _{n=0}^{{N}_{\mathrm{eq}}}{\left(\tilde{C}({t}_{n})-f({t}_{n})\right)}^{2}}.\end{eqnarray} \tag{ 11 }$$

Here, the respective functions are evaluated at available points in time t_n = n δ t, where δ t is the time step used to solve the equation of motion. The upper bound N_eq corresponds to a time at which the dynamics in question has seemingly equilibrated [19].

To get a feeling of which numerical value of constitutes a good or bad fit, we refer to the exemplary numerical data displayed in figures 3, 4, 8 and 9.

We introduce a second quantifier σ that measures how strongly the unperturbed dynamics are altered due to the perturbation in the first place.

$$\begin{eqnarray}&&\sigma =\sqrt{\displaystyle \frac{1}{{N}_{\mathrm{eq}}}\sum _{n=0}^{{N}_{\mathrm{eq}}}{\left(\tilde{C}({t}_{n})-C({t}_{n})\right)}^{2}}.\end{eqnarray} \tag{ 12 }$$

The construction of this quantity is similar to the construction of the quantity that measures the fit quality.

In the first round of our stability investigation, an exponential decay competes against a Gaussian decay. Following the strategy laid out in the previous sections, the first step is to find suitable Lanczos coefficients that comply with the four conditions presented in section 3.2.

We begin with the Gaussian decay. The Lanczos coefficients that exactly correspond to Gaussian decay of the form $C(t)={e}^{-{t}^{2}/2}$ are analytically known [method one in section 3.1)], i.e., ${b}_{n}=\sqrt{n}$ [20]. Nevertheless, the second condition, the compliance with the operator growth hypothesis, needs to be satisfied. To this end, we choose a cutoff-point n^⋆ at which the coefficients are continued in a linear fashion. This yields the following Lanczos coefficients ['g' = 'Gaussian']

$$\begin{eqnarray}{b}_{n}^{{\rm{g}}}=\left\{\begin{array}{ll}\sqrt{n}, & n\leqslant {n}^{\star }\\ \alpha n+\gamma , & n\gt {n}^{\star }.\end{array}\right.\end{eqnarray} \tag{ 13 }$$

The parameters α and γ are determined by the requirement of a smooth transition from square-root growth to linear growth, i.e., $\alpha =1/(2\sqrt{{n}^{\star }})$ and $\gamma =\sqrt{{n}^{\star }}/2$. The change from square-root to linear growth for n > n^⋆ does not strongly affect the 'Gaussianity' of the correlation dynamics [compare figure 2] such that condition one remains fulfilled.

Next, coefficients for an (approximate) exponential decay need to be found. As mentioned, a correlation function can never be truly exponential, since it necessarily features zero slope at t = 0. Thus, we only require the exponential decay to be present after a short Zeno-time, which is usually exceedingly short compared to the relaxation time [21, 22]. We achieve an approximate exponential decay via an educated guess [method two in section 3.1]. Slow dynamics are characterized by a relatively small first coefficient b₁, followed by a jump to a larger remainder of coefficients b_n≥2. Thus, we make the following ansatz for the coefficients ['e' = 'exponential']

$$\begin{eqnarray}{b}_{n}^{{\rm{e}}}=\left\{\begin{array}{ll}a, & n=1\\ \alpha n+\gamma , & n\geqslant 2.\end{array}\right.\end{eqnarray} \tag{ 14 }$$

The parameter a is set to 1.2, a similar magnitude as the first Gaussian coefficient ${b}_{1}^{{\rm{g}}}$. The other parameters α and γ are the same as in the Gaussian case. Thus, condition two, i.e., compliance with the universal operator growth hypothesis, is fulfilled. The Lanczos coefficients as defined in equation (13) and equation (14) are depicted in figure 1. Since both sets of coefficients coincide for n > n^⋆, condition four is satisfied as the quantity $q({b}_{n}^{{\rm{g}}},{b}_{n}^{{\rm{e}}})=0.999993$ is sufficiently close to unity. The corresponding dynamics are depicted in figure 2. The ansatz for the coefficients ${b}_{n}^{{\rm{e}}}$ in equation (14) indeed yields a nice exponential decay. A fit of the form Ae^{−λ t} with A = 1.02 and λ = 0.24 captures the dynamics quite well. Thus, condition one is sufficiently fulfilled for the exponential decay as well.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Further, we can extract from figure 2 that both dynamics decay on more or less the same time scale. Therefore, we also view condition three as satisfied. If any, the Gaussian dynamics is faster and thus less prone to perturbations. Now that we have checked all four conditions presented in section 3.2, we are ready to apply the perturbation as laid out in section 3.3. We set λ = 0.5. Three exemplary perturbed dynamics with respective fits for the exponential case are depicted in figure 3.

Zoom In Zoom Out Reset image size

It is evident that the perturbed dynamics are still nicely described by an exponential, only the decay constant changes. For the Gaussian decay, three exemplary perturbed dynamics with respective fits are depicted in figure 4. Two displayed perturbed Gaussian curves feature oscillations, which can not possibly be captured by a Gaussian fit ansatz. The three depicted fits are much worse than in the exponential case. Insets show corresponding Lanczos coefficients.

Zoom In Zoom Out Reset image size

The impression from these six exemplary curves is confirmed in figure 5, which shows histograms of deviations of the fits from the respective perturbed dynamics. The symbol Ω denotes the number of values within a bin of size 5 · 10⁻⁴. There is a clear division between the exponential cases (blue) and the Gaussian cases (red).

Zoom In Zoom Out Reset image size

The deviations of the exponential decays are much smaller than those of the Gaussian decays. In particular, the mean value ${\bar{\epsilon }}_{{\rm{e}}}=0.002$ (indicated by a blue, dashed line) in the exponential case is about twenty times smaller than in the Gaussian case ${\bar{\epsilon }}_{{\rm{g}}}=0.042$. For a visualization of this comparison see figures 3 and 4, whose exemplary curves feature values of close to the respective mean values.

We have to be aware that this apparent stability of the exponential may be due to the possibility that exponential decays are generally less affected by our constructed perturbation than the Gaussian decays. To check this possibility and to show that this is not the case, we consider the quantifier σ, which measures the difference between the perturbed and the unperturbed dynamics. The inset of figure 5 depicts a scatter plot of all 1000 pairs (σ_i , _i ). If one decay would be consistently less affected than the other, a cluster of points to the left [small values of σ] would emerge. This is evidently not the case, as the distribution along the horizontal axis is relatively similar for both decays. The existence of the diagonal edge is expected, since there can be no fit that is 'farther away' from the perturbed dynamics than the unperturbed dynamics itself, which is, of course, also an eligible candidate for fitting.

Thus, we conclude that exponential decays is indeed stable with respect to perturbations [as in equation (10)]. In contrast, Gaussian decays seem to be quite unstable. This is the first main result of the paper at hand.

In the second round, we investigate and compare two types of damped oscillations, one with an exponential damping factor, the other with a Gaussian damping factor. As in the last section, first, two sets of Lanczos coefficients need to be found that satisfy all four conditions imposed in section 3.2.

We again begin with the Gaussian case. The coefficients ${b}_{n}^{\mathrm{gdo}}$ ['gdo' = 'Gaussian damped oscillation'] that correspond to an oscillation damped by a Gaussian factor a neither analytically known, nor can we make an educated guess. This only leaves the third method, in which we 'reverse-engineer' the coefficients from the dynamics itself. To this end, we choose a particular correlation function $C(t)={e}^{-{t}^{2}/8}\cos (2t)$ and proceed as detailed in appendix A. This procedure allows us to obtain about 50 Lanczos coefficients. After that, the coefficients are continued 'by hand' in a manner that first, respects the pattern exhibited by the coefficients so far, and second, eventually becomes linear.

For the exponentially damped oscillation we can again make an educated guess [method two] for the coefficients ${b}_{n}^{\mathrm{edo}}$ ['edo' = 'exponentially damped oscillation'].

We set the first two coefficients b₁, b₂ to some values (b₁ = 2.0, b₂ = 1.6), followed by a jump to larger coefficients b_n≥3, which then grow in a linear fashion. The slope is determined by and coincides with the slope in the Gaussian case.

The unperturbed Lanczos coefficients and corresponding correlation functions are depicted in figure 6 and figure 7, respectively. The 'reverse-engineered and continued by hand' coefficients ${b}_{n}^{\mathrm{gdo}}$ indeed reproduce the given correlation function $C(t)={e}^{-{t}^{2}/8}\cos (2t)$. Further, the guess for ${b}_{n}^{\mathrm{edo}}$ yields a nice exponentially damped curve, i.e., the fit ansatz $\tilde{C}(t)={{Ae}}^{-\mu t}\cos (\omega t-\phi )$ with A = 1.04, μ = 0.57, ω = 2.19 and ϕ = − 0.32 captures the dynamics quite well. Therefore, condition one is satisfied for both cases. The two sets of coefficients were designed to eventually attain linear growth, thus condition two is fulfilled. It is evident from figure 7 that both dynamics decay on similar time scales, rendering condition three satisfied. Lastly, the quantity $q({b}_{n}^{\mathrm{gdo}},{b}_{n}^{\mathrm{edo}})$ practically equals unity. Thus, all four conditions are met.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Next, we switch on the perturbation and set λ = 0.1. Exemplary perturbed dynamics are displayed in figures 8 and 9, respectively. The perturbed dynamics are fitted with the 4-parametric ansatz $\tilde{C}(t)={{Ae}}^{-\mu t}\cos (\omega t-\phi )$ in the exponential and $\tilde{C}(t)={{Ae}}^{-\mu {t}^{2}}\cos (\omega t-\phi )$ in the Gaussian case. While the perturbed coefficients in the exponential case still lead to correlation functions that are within the class of exponentially damped oscillations [with different decay constants μ and frequencies ω], the same can not be said for the Gaussian case, where the perturbed dynamics decay too slowly and can not by captured by the fit ansatz above.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The impression from the six exemplary curves is confirmed in figure 10, which displays histograms of the quantifier . Again, the division into two distinct peaks is evident. The respective mean values ${\overline{\epsilon }}_{\mathrm{gdo}}=0.026$ and ${\overline{\epsilon }}_{\mathrm{edo}}=0.006$ are indicated as dashed, vertical lines and differ by a factor of about five. For a visualization of what these values mean for the fit quality, see figures 8 and 9, whose curves feature values of close to the respective mean values. Thus, we conclude that the exponentially damped dynamics are quite stable and, in particular, more stable than the oscillations damped by a Gaussian. The disparity between the two is not as strong as for the decays in the previous section [however, different values of λ may not be directly comparable.]

Zoom In Zoom Out Reset image size

The inset of figure 10 again shows a scatter plot of all points (σ_i , _i ). The blue cluster of dots is a little more concentrated to small values of σ, indicating that the exponentially damped oscillations are a little less altered by the perturbation than their Gaussian counterpart. However, the marginal distribution over σ is still less partitioned into two peaks than the one over [which is the data in the histogram itself]. Hence, the stability of exponentially damped oscillations is still due to the nature of the dynamics and not due to a smaller effect of the perturbation on the dynamics.

We conclude that exponentially damped oscillations are stable with respect to perturbations [as in equation (10)]. In contrast, oscillations damped by a Gaussian seem to be quite unstable. This is the second main result of the paper at hand.

In this section, we lift the restriction imposed on the v_n earlier, which was that frequencies were cut at N_f ≈ d/3. Instead, we include all frequencies in the construction of the perturbation, which amounts to setting N_f = d. This case simply corresponds to adding uncorrelated random numbers to the unperturbed coefficients b_n . Otherwise, the strategy is pursued in the same manner as above.

Three exemplary dynamics for exponential decays are depicted in figure 11 [for conciseness, we refrain from showing more exemplary data for the other cases]. As is evident, the exponential decay is completely absent and replaced by quite irregular dynamics, which do not seem to reach an equilibrium [N_eq is set as large as possible [19] in these cases]. These curves hint at the presence of localization effects, which prevent the 'particle' [in the tight-binding picture] to leave the first site.

Zoom In Zoom Out Reset image size

A histogram of the Gaussian versus exponential data can be viewed in figure 12. The mean values read ${\bar{\epsilon }}_{{\rm{g}}}=0.15$ in the Gaussian case and ${\bar{\epsilon }}_{{\rm{e}}}=0.12$ in the exponential case, which is some orders of magnitude larger than before when shorter wavelengths were excluded. Further, the accumulation of points along the diagonal in the inset suggests that the curves no longer resemble their original form, which is expected judging from figure 11.

Zoom In Zoom Out Reset image size

The histogram for the oscillating cases is displayed in figure 13. The mean deviations read ${\bar{\epsilon }}_{\mathrm{gdo}}=0.035$ and ${\bar{\epsilon }}_{\mathrm{edo}}=0.025$. These values are comparable to the Gaussian case in figure 10. However, since the quantifier measures more or less the deviation 'per time step', the comparison of values $\bar{\epsilon }$ between dynamics that do and do not equilibrate can be void of meaning. In practice, judging from all figures displaying exemplary curves, a value of about one percent corresponds to fits that 'look good to the eye'.

Zoom In Zoom Out Reset image size

Based on these observations, we specify a property of perturbations that lead to non-generic or pathological dynamics. Namely, a perturbation is 'untypical', if it gives rise to Lanczos coefficients that vary non-smoothly, i.e., there is no minimal correlation length within the coefficients. Again, what this means on the level of Hermitian matrices is difficult to assess and beyond the scope of this work. This is the third main result of the paper at hand.